A personalized and non-parametric framework for detecting changes in gait cycles

Department of Mechanical & Aerospace Engineering, University at Buffalo

Department of Mechanical & Aerospace Engineering, University at Buffalo

Department of Industrial and Systems Engineering, University at Buffalo

Farmer School of Business, Miami University

Farmer School of Business, Miami University

---

1 R Setup and Required Packages

In the following code chunk we load the packages used to support our analysis.

# check if packages are not installed; if yes, install missing packages
packages = c("tidyverse", "magrittr", # typical data analysis packages
             "MALDIquant", # match closest points between two vectors
             "foreach", "doParallel", "parallel", # packages for parallelization
             "R.matlab")
newPackages = packages[!(packages %in% installed.packages()[,"Package"])]
if(length(newPackages) > 0) install.packages(newPackages)

# using the library command to load all packages; invisible used to avoid printing all packages and dependencies used
invisible(lapply(packages, library, character.only = TRUE))

# source("./Functions.R") # our custom built functions

# set.seed(2020)
# startTime <- Sys.time()

---

2 Gait Acceleration Data and Rational Subgroups

We converted the raw IMU acceleration signals in the local to the global reference frame and removed the gravity. The acceleration signals were then transformed back to the local reference frame. Sagittal acceleration, lateral acceleration, and acceleration magnitude signals were then calculated. We used the vertical acceleration component to segment the gait cycles in order to isolate individual gait cycles. The acceleration profiles along with the experimental time stamps were stored in mat files.

In this code chunk we load the segmented acceleration mat files. The \(1^{st}\) 10 minutes of the data were considered as warm up period and thus excluded form the analysis. The acceleration profiles of the gait cycles during each walking cycle were grouped into a rational subgroup. We stored the start and end times of the walking cycles in csv files and loaded and used them in the following code chunk. The acceleration profiles within each start and end times to subgroup the gait cycles.

for (id in setdiff(1:15, 13)) {
  #################### Read start and end time of the subgroups
  video <- read.csv(paste0(file="../Data/csvFiles/Sub", id, ".csv"))

  video_leave <- video$Leaves
  video_leave <- video_leave[!is.na(video_leave)]
  video_leave <- video_leave - 600
  video_leave <- video_leave[video_leave > 0]

  video_enter <- video$Enters
  video_enter <- video_enter[!is.na(video_enter)]
  video_enter <- video_enter - 600
  video_enter <- video_enter[video_enter > 0]

  ################################################### Read from mat files
  raw1 <- readMat(paste0("../Data/matFiles/Subject", id, "_aZ_seg.mat"))

  raw2 <- raw1$gait

  num_rows <- length(raw2)/5

  aM <- list("vector")
  aS <- list("vector")
  aL <- list("vector")
  exp_time <- c()
  for (i in 1:num_rows) {
    aM[[i]] <- raw2[[1*num_rows + i]][[1]][1,]
    aS[[i]] <- raw2[[2*num_rows + i]][[1]][1,]
    aL[[i]] <- raw2[[3*num_rows + i]][[1]][1,]
    exp_time[i] <- raw2[[i+4*num_rows]][[1]][1,1]
  }

  ###################################################### Create subgroups
  if (video_leave[1] > video_enter[1]){
    video_leave <- c(exp_time[1], video_leave)
  }

  if (tail(video_leave, 1) > tail(video_enter, 1)){
    video_enter <- c(video_enter, tail(exp_time, 1))
  }

  conf <- 0.05 * mean(video_enter - video_leave) # 5% confidence for each subgroup

  # Match video subroup times against the experimental times from mat files
  leave_index <- match.closest(video_leave + conf, exp_time)
  enter_index <- match.closest(video_enter - conf, exp_time)

  aM_list <- list()
  aS_list <- list()
  aL_list <- list()
  t_list <- list()
  for (i in 1:length(leave_index)) {
    aM_list[[i]] <- aM[leave_index[i]:enter_index[i]]
    aS_list[[i]] <- aS[leave_index[i]:enter_index[i]]
    aL_list[[i]] <- aL[leave_index[i]:enter_index[i]]
    t_list[[i]] <- exp_time[leave_index[i]:enter_index[i]]
  }

  len <- c()
  for (i in 1:length(aM_list)) {
    len <- c(len, lengths(aM_list[[i]]))
  }

  cut_len <- quantile(len, 0.02)

  for (i in 1:length(aM_list)) {
    remove_ind <- which(lengths(aM_list[[i]]) < cut_len)
    remove_ind <- c(remove_ind, 3000)
    aM_list[[i]] <- aM_list[[i]][-remove_ind]
    aS_list[[i]] <- aS_list[[i]][-remove_ind]
    aL_list[[i]] <- aL_list[[i]][-remove_ind]
    t_list[[i]] <- t_list[[i]][-remove_ind]
  }

  batchSize <- lengths(aM_list)

  ###############################################
  aM_list_all <- do.call(c, aM_list)
  aS_list_all <- do.call(c, aS_list)
  aL_list_all <- do.call(c, aL_list)
  t_all <- do.call(c, t_list)
  ###############################################
  aM_mat <- do.call(rbind, aM_list_all)[,1:cut_len]
  aS_mat <- do.call(rbind, aS_list_all)[,1:cut_len]
  aL_mat <- do.call(rbind, aL_list_all)[,1:cut_len]

  ##################
  assign(paste0("sub", id, "_aMSLT"),
         list(aMag=aM_mat, aSag=aS_mat, aLat=aL_mat, time=t_all, batchSize=batchSize))
}

########################################### Save
save(sub1_aMSLT, sub2_aMSLT, sub3_aMSLT, sub4_aMSLT, sub5_aMSLT, sub6_aMSLT,
     sub7_aMSLT, sub8_aMSLT, sub9_aMSLT, sub10_aMSLT, sub11_aMSLT, sub12_aMSLT,
     sub14_aMSLT, sub15_aMSLT,

     file="../Data/rData/subGs_for_depth.Rdata")

3 Depth Calculation

Explain about depth concept in general; centrality and outlyingness

3.1 Mode Depth

Calculate mode depth…

load(file="../Data/rData/subGs_for_depth.Rdata")

cores = detectCores() - 2 # to give the server some breathing Room
cl = makePSOCKcluster(cores)
registerDoParallel(cl)

for (id in setdiff(1:15, 13)) {
  
  aMag <- get(paste0("sub", id, "_aMSLT"))$aMag
  tExp <- get(paste0("sub", id, "_aMSLT"))$time
  
  ################################## Mode Depth
  inControl <- aMag[1:500,]
  inControlDepth <- {depth.mode(fdata(inControl))}$dep
  
  ###################### foreach
  onlineDepth <- c()
  end <- nrow(aMag)
  onlineDepth <- foreach(i=501:end, .packages = c('fda.usc', 'tidyverse'),
                         .combine='c') %dopar% {
                           temp1 <- aMag %>% .[i,]
                           append1 <- rbind(inControl, temp1)
                           temp2 <- depth.mode(fdata(append1))
                           
                           temp3 <- temp2$dep %>% .[501]
                           temp3
                         }
  
  modeDepth <- c(inControlDepth, onlineDepth)
  
  assign(paste0("sub", id, "_mode_mag"),
         list(modeDepth=modeDepth, tExp=tExp))
}

save(sub1_mode_mag, sub2_mode_mag, sub3_mode_mag, sub4_mode_mag, sub5_mode_mag, sub6_mode_mag,
     sub7_mode_mag, sub8_mode_mag, sub9_mode_mag, sub10_mode_mag, sub11_mode_mag, sub12_mode_mag,
     sub14_mode_mag, sub15_mode_mag,
     file="../Data/rData/mode_mag.Rdata")

3.2 Visualize the Mode Depth

load(file="../Data/rData/mode_mag.Rdata")

for (id in setdiff(1:15, 13)) {
  depthMode <- get(paste0("sub", id, "_mode_mag"))$modeDepth
  tMode <- get(paste0("sub", id, "_mode_mag"))$tExp

  cat("###", paste0("Subject", id, " Mode Depth"), "{-}",'\n')
  plot(tMode, depthMode, pch=16, cex=1,col=c(rep("red", 500), rep("black", length(depthMode)-500)))
  legend("topright", legend=c("Baseline data", "New data"), col=c("red", "black"))
  cat('\n \n')

}

Subject1 Mode Depth

Subject2 Mode Depth

Subject3 Mode Depth

Subject4 Mode Depth

Subject5 Mode Depth

Subject6 Mode Depth

Subject7 Mode Depth

Subject8 Mode Depth

Subject9 Mode Depth

Subject10 Mode Depth

Subject11 Mode Depth

Subject12 Mode Depth

Subject14 Mode Depth

Subject15 Mode Depth

3.3 MFHD Depth

Calculate MFHD depth…

Explain since it take very long it was performed on supercomputer cluster in parallel and the eval option in the following code chunk was set to false since this could take weeks to run on the local computer.

#################################################
load(file="../Data/rData/subGs_for_depth.Rdata")

cores = detectCores() - 2 # to give the server some breathing Room
cl = makePSOCKcluster(cores)
registerDoParallel(cl)

for (id in setdiff(1:15, 13)) {
  
  aSag <- get(paste0("sub", id, "_aMSLT"))$aSag
  aLat <- get(paste0("sub", id, "_aMSLT"))$aLat
  tExp <- get(paste0("sub", id, "_aMSLT"))$time
  
  ################################## Mode Depth
  inControlSag <- aSag[1:500,]
  inControlLat <- aLat[1:500,]
  inControlDepth <- {MFHD(y1=inControlSag, y2=inControlLat, alpha=0.125, Beta=0.5)}$MFHDdepth[1,]
  
  ###################### foreach
  onlineDepth <- c()
  end <- nrow(aSag)
  onlineDepth <- foreach(i=501:end, .packages = c('fda.usc', 'tidyverse', 'MFHD'),
                         .combine='c') %dopar% {
                           temp1 <- aSag %>% .[i,]
                           temp2 <- aLat %>% .[i,]
                           
                           append1 <- rbind(inControlSag, temp1)
                           append2 <- rbind(inControlLat, temp2)
                           
                           temp3 <- MFHD(y1=append1, y2=append2, alpha=0.125, Beta=0.5)
                           
                           temp4 <- temp3$MFHDdepth %>% .[1,501]
                           temp4
                         }
  
  MFHDdepth <- c(inControlDepth, onlineDepth)
  
  assign(paste0("sub", id, "_MFHD"),
         list(MFHDdepth=MFHDdepth, tExp=tExp))
}

save(sub1_MFHD, sub2_MFHD, sub3_MFHD, sub4_MFHD, sub5_MFHD, sub6_MFHD, sub7_MFHD, sub8_MFHD,
     sub9_MFHD, sub10_MFHD, sub11_MFHD, sub12_MFHD, sub14_MFHD, sub15_MFHD,
     file="../Data/rData/MFHD.Rdata")

3.4 Visualize the MFHD Depth

load(file="../Data/rData/MFHD.Rdata")

for (id in setdiff(1:15, 13)) {
  depthMFHD <- get(paste0("sub", id, "_MFHD"))$MFHDdepth
  tMFHD <- get(paste0("sub", id, "_MFHD"))$tExp
  
  cat("###", paste0("Subject", id, " MFHD Depth"), "{-}",'\n')
  plot(tMFHD, depthMFHD, pch=16, cex=1,col=c(rep("red", 500), rep("black", length(depthMFHD)-500)))
  legend("topright", legend=c("Baseline data", "New data"), col=c("red", "black"),
         pch=16)
  cat('\n \n')
  
}

Subject1 MFHD Depth

Subject2 MFHD Depth

Subject3 MFHD Depth

Subject4 MFHD Depth

Subject5 MFHD Depth

Subject6 MFHD Depth

Subject7 MFHD Depth

Subject8 MFHD Depth

Subject9 MFHD Depth

Subject10 MFHD Depth

Subject11 MFHD Depth

Subject12 MFHD Depth

Subject14 MFHD Depth

Subject15 MFHD Depth

---

References

---

1. Email: saebraga@buffalo.edu | Website: LinkedIn↩︎

2. Email: jiyeonk@buffalo.edu | Phone: +1-716-645-6063 | Website: University at Buffalo Official↩︎

3. Email: loracavu@buffalo.edu | Phone: +1-716-645-4696 | Website: University at Buffalo Official↩︎

4. Email: fmegahed@miamioh.edu | Phone: +1-513-529-4185 | Website: Miami University Official↩︎

5. Email: farmerl2@miamioh.edu | Phone: +1-513-529-4823 | Website: Miami University Official↩︎